# 2013/uns13/resources: optsparse.py

File optsparse.py, 3.6 KB (added by drasmuss, 5 years ago) |
---|

Line | |
---|---|

1 | """ |

2 | Implementation of sparse connectivity for NEF/SpiNNaker |

3 | The Spinn constructor can be passed a max_fan_in parameters, |

4 | which limits the number of connections a neuron can receive in the network. |

5 | |

6 | Author: Terry Stewart |

7 | """ |

8 | |

9 | import numeric as np |

10 | import random |

11 | import ca.nengo.util.MU as MU |

12 | #import Jama.Matrix as Matrix |

13 | import subprocess |

14 | import os |

15 | import struct |

16 | from array import array |

17 | from java.lang import System |

18 | |

19 | #uncomment this function to calculate the pseudoinverse directly within |

20 | #nengo's jython |

21 | #def pinv(matrix, limit): |

22 | # m = Matrix([[x for x in row] for row in matrix]) |

23 | # svd = m.svd() |

24 | # sInv = svd.getS().inverse() |

25 | # |

26 | # i = 0 |

27 | # while i<svd.getS().rowDimension and svd.getS().get(i, i)>limit: |

28 | # i += 1 |

29 | # while i<len(matrix): |

30 | # sInv.set(i, i, 0) |

31 | # i+=1 |

32 | # result = svd.getV().times(sInv).times(svd.getU().transpose()).getArray() |

33 | # return result |

34 | |

35 | #this function calculates the pseudoinverse through an external callout |

36 | #to a python file, allowing the use of numpy |

37 | def pinv(matrix, limit): |

38 | shortfile = "matrix_%f" % random.random() |

39 | filename = os.path.join("external", shortfile) |

40 | |

41 | #write matrix data to file |

42 | infile = file(filename, "wb") |

43 | float_array = array('f', [x for row in matrix for x in row]) |

44 | float_array.tofile(infile) |

45 | infile.close() |

46 | |

47 | #create output file |

48 | outfile = file(filename+".inv", "wb") |

49 | outfile.close() |

50 | |

51 | # print System.getProperty("os.name") |

52 | if "windows" in System.getProperty("os.name").lower(): |

53 | subprocess.call("cmd /c pseudoInverse.bat "+shortfile+" "+shortfile+".inv"+" "+str(limit)+" "+str(-1), |

54 | stdout=subprocess.PIPE, stderr=subprocess.PIPE, cwd="external") |

55 | else: |

56 | subprocess.call("external/pseudoInverse "+filename+" "+filename+".inv"+" "+str(limit)+" "+str(-1), shell=True) |

57 | |

58 | mfile = file(filename+".inv", "rb") |

59 | invmatrix = array("f") |

60 | invmatrix.fromfile(mfile, len(matrix)*len(matrix[0])) |

61 | return [[invmatrix[i*len(matrix[0])+j] for j in range(len(matrix[0]))] for i in range(len(matrix))] |

62 | |

63 | def compute_sparse_weights(origin, post, transform, fan_in, noise=0.1, num_samples=100): |

64 | encoder = post.encoders |

65 | radius = post.radii[0] |

66 | |

67 | if hasattr(transform, 'tolist'): transform=transform.tolist() |

68 | |

69 | approx = origin.node.getDecodingApproximator('AXON') |

70 | |

71 | # create X matrix |

72 | X = approx.evalPoints |

73 | X = MU.transpose([f.multiMap(X) for f in origin.functions]) |

74 | |

75 | # create A matrix |

76 | A = approx.values |

77 | |

78 | S = fan_in |

79 | N_A = len(A) |

80 | samples = len(A[0]) |

81 | N_B = len(encoder) |

82 | w_sparse = np.zeros((N_B, N_A),'f') |

83 | noise_sd = MU.max(A)*noise |

84 | num_samples = min(num_samples,N_B) |

85 | decoder_list = [None for _ in range(num_samples)] |

86 | for i in range(num_samples): |

87 | indices = random.sample(range(N_A), S) |

88 | activity = [A[j] for j in indices] |

89 | n = [[random.gauss(0, noise_sd) for _ in range(samples)] for j in range(S)] |

90 | activity = MU.sum(activity, n) |

91 | activityT = MU.transpose(activity) |

92 | gamma = MU.prod(activity, activityT) |

93 | |

94 | upsilon = MU.prod(activity, X) |

95 | |

96 | gamma_inv = pinv(gamma, noise_sd*noise_sd) |

97 | |

98 | |

99 | decoder_list[i] = MU.prod([[x for x in row] for row in gamma_inv], upsilon) |

100 | |

101 | for i in range(N_B): |

102 | ww = MU.prod(random.choice(decoder_list), MU.prod(MU.transpose(transform), encoder[i])) |

103 | |

104 | for j, k in enumerate(indices): |

105 | w_sparse[i,k] = float(ww[j])/radius |

106 | |

107 | return list(w_sparse) |

108 |