text:
# Reduced rank regression

#1

Hello Stan Community,

I have been attempting to fit a reduced-rank regression model with Stan. Briefly, reduced-rank regression takes the multivariate regression Y=XB+e with p predictors (in X) and q outcomes (in Y), and decomposes the p-by-q coefficient matrix B as B = UV', with U being p-by-k and V being q-by-k, where k < min(p, q).

This general form introduces a couple of complications, in that different U and V can result in identical likelihoods. More specifically, flipping the sign on both U and V and permuting the columns of U and V produce equivalent models.

I have managed to get a model where k=1 to work to my satisfaction [it recovers parameters of simulated data and produces interpretable results on real data] by using a transformed parameter that forces the signed sum of squares of V to be positive. It would be simpler to constrain just one entry of V to be positive, but I have no good reason to expect a particular value of V to be positive, and with the data I have so far, none of the posterior distributions of the entries in V are far enough from zero to not fold some non-trivial part of the tail over the zero line.

functions {
int sign(real x) {
if (x == 0)
return 0;
else if (x < 0)
return -1;
else
return 1;
}
}
data {
int N; //number of rows
int P; //number of predictors
int C; //number of outcomes
real<lower=0> sig; //sigma of prior t
real<lower=1> nuf; //nu of priors on coefficients
matrix[N, C] y; //outcome
matrix[N, P] X; //design matrix
}
parameters {
real<lower=1> nu; //nu of noise student t distribution
vector[P] urU;
vector[C] urV;
vector<lower=0>[C] sigma;
}
transformed parameters {
vector[P] U;
vector[C] V;
{ //scope hider
//point V in a positive direction
real s;
real sgn;
s = 0;
for (i in 1:C) {
s += sign(urV[i])*urV[i]^2;
}
sgn = sign(s);
for (i in 1:C)
V[i] = urV[i] * sgn;
for (j in 1:P)
U[j] = urU[j] * sgn;
} // end scope hider
}
model {
nu ~ gamma(2, 0.1);
sigma ~ cauchy(0, 1);
urU ~ student_t(nuf, 0, sig);
urV ~ student_t(nuf, 0, sig);
{
matrix[N,C] mu;
mu = X*U*V';
for (j in 1:C)
y[:,j] ~ student_t(nu, mu[:,j], sigma[j]);
}
}

Anyway, things really fall apart when k>1. I suspect it is because of the column-ordering issue, but there may be some additional problem that I am missing. I have tried introducing a k by k diagonal matrix S (i.e. Y=XUSV') with the diagonal entries constrained to be positive and ordered, but that hasn’t really helped. A while ago I looked at orthogonalizing the columns of U with the Gram-Schmidt process, but I have learned a lot about Stan since then, so I should perhaps revisit that idea. Here is the Stan model for k>1 (that doesn’t really work and doesn’t use the ordered S). Apologies for changing variable names, here Y=XBA'+e, and B and A have qstar columns.

functions {
int sign(real x) {
if (x == 0)
return 0;
else if (x < 0)
return -1;
else
return 1;
}
}
data {
int N; //number of rows
int P; //number of predictors
int<lower=2> C; //number of outcomes
int<upper=min(C,P)> qstar; //ceiling on rank
matrix[N, C] y; //outcome
matrix[N, P] X; //design matrix
}
parameters {
matrix[P, qstar] B0;
matrix[C, qstar] A0;
vector<lower=0>[C] sigma;
matrix<lower=0>[P, qstar] lam;
positive_ordered[qstar] tau;
real<lower=0> nu;
}
transformed parameters {
matrix[P, qstar] B;
matrix[C, qstar] A;
for (j in 1:P) {
for (h in 1:qstar) {
B[j,h] = B0[j,h]*lam[j,h]*tau[h];
}
}
A = A0;
//point A in a generally positive direction
for (c in 1:qstar) {
real ss = 0;
for (r in 1:C)
ss += sign(A[r, c])*A[r, c]*A[r, c];
A[:,c] = A[:,c]*sign(ss);
B[:,c] = B[:,c]*sign(ss);
}
}
model {
matrix[N, C] mu;
nu ~ gamma(2, 0.1);
sigma ~ cauchy(0, 1); // different from reference paper
to_vector(A0) ~ normal(0, 1);
to_vector(B0) ~ normal(0, 1);
to_vector(lam) ~ cauchy(0, 1);
tau ~ cauchy(0, 1);
mu = X*B*A';
for (i in 1:C)
y[:,i] ~ student_t(nu, mu[:,i], sigma[i]);
}

As a rank newbie to Stan (and Bayesian modelling in general), I would appreciate any help figuring out how to get something like rank-reduced regression working. I know there are a lot of similar approaches with different names, but in the papers I have found describing related methods, I get lost before I understand how the sign and column permutation problems get solved.

I’d also welcome any general tips or advice on my models. Thanks!

#2

We don’t usually recommend this other than in doing something like a QR decomposition of the predictors.

If you do want to do this, then I would worry about correcting for sign after the fact. You probably only care about preiction anyway with a technique like this.

You can simplify:

real s = sign(urV) .* square(urV);
vector[C] V = urV * sgn;
vector[C] U = urU * sgn;

If what you want to do is a kind of variable selection, then I’d recommend the recent work of Piironen and Vehtari on horsehoe priors.

#3

Reduced rank regression is not something you see Bayesians do very often, but in order to get around problems you mention you can use a rectangular Cholesky factor for \mathbf{U}, as in

data {
cholesky_factor_cov[p, k] U;
...

The columns cannot be swapped without destroying the lower-trapezoidal structure and the signs cannot be flipped without making the diagonal elements negative.

#4

Thanks for the thoughtful replies!

So it sounds like my approach is pretty unusual, but I suspect that my problem is not that unusual, so maybe there is a more conventional approach to the problem:

I have a dataset from an experiment with about 10 outcome measures per case and a number of possible predictors/covariates that approaches the number of cases. The expectation is that there will be correlations among the predictors and among the outcomes, and I expect that some groups of the correlated predictors will be related to some groups of correlated outcomes. My goal is to identify relationships between predictors and outcomes so I might use those relationships in future studies and practical applications.

For example, I want to be able to say that if you care about outcome measure 4, then you should find people who score high on these three predictors and low on these other two, and by the way, you can then expect for outcomes 5-7 to be relatively high, but 8 will probably be low.

So the question is whether there exists a conventional approach to solve problems like this in a Bayesian framework.

Their papers, along with the related case study of Michael Betancourt, Bayes Sparse Regression have been helpful. The examples (unless I’m missing something) seem focused on univariate models, and it isn’t clear to me how to extend them to multivariate cases. If I assume there is one pattern or group of predictors that accounts for one pattern or group of outcomes, then the results I get make sense, but I suspect I might be missing possible explanations for the other outcomes that don’t fall in that pattern.

If I understand this suggestion, I would be constraining the model such that the ordering of the columns of X now becomes an important decision. The second and subsequent latent factors in X would be prohibited from having any contribution from column 1, because of the zeros in the upper triangular part of U. I might be able to come up with a reason to justify some particular column order for X. I’ll certainly give it a try.

It might be the case that there is no perfect solution, here. Thanks again for the suggestions so far (and the vectorization tip).

#5

#6

For rows 2 … k, there are exclusion restrictions in the previous columns, so yes the order is relevant. But I don’t get the sense that people who are doing reduced rank regression are particularly interested in interpreting particular coefficients. The elements of V\prime can adjust to still fit the first k variables well.

#7

But it seems that’s not the case?

Our projection predictive variable selection code does not currently take into account the correlations, although if used with latent variable models it might be usable for making some latent variables to shrink near zero.

Let’s be careful with word multivariate as many would say that many covariates and one target is also multivariate model. We demonstrate only with univariate targets, but it’s easy to extend to multivariate targets, but it doesn’t matter if you are more interested in the reduced dimensional latent linear model than variable selection.

So it seems I can’t help with the reduced rank part, and it seems you don’t need yet the variable selection part.

#8

Reframing your problem, you are essentially trying to infer a factor model of your latent regression weight matrix. I spent this spring trying to get factor models in Stan to work. Due to the non-identifiability UV^T=U(R^{-1}R)V^T, you need to impose a lot of constraints, as you have already started doing. I would say that you can get it to work if all you want is a low-rank decomposition i.e. some arbitrary, but consistent, pair \hat{U} and \hat{V}, but I would say that it is almost impossible to get e.g. PCA vectors with uncertainty quantification unless their effective singular values (your S) are very well-separated.

The gold standard is probably triangular constraints like used in this blog post by Rick Farouni. This article details how you can construct a proper prior for that setting, equivalent to the standard Gaussian choice.

If you want to try something with orthogonal matrices, I can dig around for the forum entry where we discuss a suitable design using Cayley transforms that works in Stan (Gram-Schmidt not necessary). But I did not really find them to be all that helpful personally since there is still the permutation issue.

Ironically, it appears that Stan’s strength makes it harder’’ to sample from these models since it can effortlessly jump between symmetry modes :)

#9

Thanks to the helpful suggestions above, I have arrived at something that appears to work on simulated data most of the time [if my simulation pulls from the tails of the cauchy priors, that sometimes louses things up]. It doesn’t seem to be anything novel, but I’m going to post my model and some comments to maybe save others (and future me) some time. If anyone spots any problems or improvements on this method, I would be very happy to hear about it. The model and some code to simulate data are at the end of the post, if you want to skip right to that.

I want interpretable reduced rank regression. Here I explain what I think that means.

To restate the goal, here: I’ve got some predictors/covariates in a matrix X, and some outcomes of interest in matrix Y. I’m interested in patterns of predictors that are related to patterns of outcomes, and I want to interpret the coefficients that describe these patterns (i.e. I’m not exclusively interested in predicting new cases). This is described in the literature as reduced rank regression, and I think partial least squares regression is a special case.

The approach is to fit a model like Y=XBA'+e where B and A are matrices with p and c rows, respectively and k columns, where p is the number of predictors and c is the number of outcomes. To try to get a sparser solution, I used a horseshoe prior for B. In order to get a stable solution to this problem with Stan, we need several constraints on the coefficients A and/or B.

I went with the gold standard. As the blog post indicates, this introduces a complication in that if you constrain B to be lower trapezoidal, then the column order of X matters, since the first column of X only contributes to the first latent factor, and the second column only to the first two, etc. And it turns out that this constraint only successfully forces the model to be unimodal if the coefficients on the diagonal are large relative to the rest of the values in the column. If the model fits about as well whether a particular coefficient is positive or negative (because its true value is close to zero), then constraining that coefficient doesn’t really constrain the model.

It also seems like it is critical to set some reasonable initial values, and that each chain needs to start with approximately the same initial values. At first, I cheated by using the simulated values as the initial values, and that works but is implausible in real data.

The answer to the column ordering and initial value problems are sort of related, in that if you can get a good initial guess, you can re-order X to give yourself a better chance at finding a stable solution.

So I came up with a way to find an approximate solution with about the right shape. Find the least-squares estimate of \beta (where Y=X\beta) and use QR decomposition to decompose \beta' to a unitary square matrix and a triangular matrix, and then take the first k rows of each to get \hat{A} and \hat{B} (transposed). Now, \hat{A} is constrained by the QR decomposition to have L2 norm = 1, which is not part of the constraint in our model. If we use the standard normal prior for A, then the expected L2 norm of A is \sqrt{c*k}, so just multiply \hat{A} by that and divide \hat{B} by the same.

Now I can inspect \hat{B} and rearrange the columns of X so that coefficients with relatively large absolute values would be moved to the diagonals of the rearranged \hat{B}. After this re-arrangement, I can do the above approximation again and confirm the new decomposition has large values on the diagonals. This second round of estimates can then serve as the initial values for the chains in Stan.

Here is the stan code I’m using. I’ve only really tested it with smallish p and c, and k<=3:

data {
int n;
int p;
int c;
int k;
real<lower=0> sigma;
matrix[n, p] X;
matrix[n, c] Y;
}
transformed data {
int ntrap = p*k - k*(k-1)/2 - k;
}
parameters {
vector<lower=0>[k] diags;
vector<lower=0>[k] diag_sig;
vector[ntrap] lowtrap;
vector<lower=0>[ntrap] lowtrap_sig;
real<lower=0> global_shrink;
matrix[c, k] A;
}
transformed parameters {
matrix[p,k] B;
{ //hide integer
int idx;
idx=0;
B = rep_matrix(0, p, k);
for (col in 1:k) {
B[col, col] = diags[col]*diag_sig[col]*global_shrink;
for (r in (col+1):p) {
idx+=1;
B[r, col] = lowtrap[idx]*lowtrap_sig[idx]*global_shrink;
}
}
} //close hider
}
model {
matrix[n, c] mu;
lowtrap_sig ~ cauchy(0,1);
diag_sig ~ cauchy(0,1);
global_shrink ~ cauchy(0,1);
lowtrap ~ normal(0,1);
diags ~ normal(0,1);
to_vector(A) ~ normal(0,1);
mu = X*B*A';
for (i in 1:c)
Y[:,i] ~ normal(mu[:,i], sigma);
}
generated quantities {
vector[n] log_lik = rep_vector(0, n);
{ //hide mu to not clutter output
matrix[n, c] mu;
mu = X*B*A';
for (row in 1:n)
log_lik[row] = normal_lpdf(Y[row,:] | mu[row,:], sigma);
}
}

And some MATLAB (2017a) code to generate data and initial values:

%% Simulate Y from known X and parameters (B, A)
p = 8; %number of columns in X
c = 6; %number of columns in Y
k = 2; %number of columns in B and A
n = 93;%number of rows in X and Y
NC = 12; %number of chains we're going to run
noise_sig = 0.1; % fixed for now.

init_sd = 0.05; % initializations use this to jitter the initial values for different chains.

% some convenient indeces
diag_idx = find(eye(p),k);
trap_idx = tril(true(p),-1);
trap_idx = trap_idx(:, 1:k);

% Select sigma B and sigma G (global)
tdist = makedist('tlocationscale', 'mu', 0, 'sigma', 1, 'nu', 1);
B_sig = abs(tdist.random(p,k));
G_sig = abs(tdist.random(1));

% give X some interesting structure so the method has something to find
sigma = vine(p,1); % a random correlation matrix, from https://stats.stackexchange.com/questions/2746/how-to-efficiently-generate-random-positive-semidefinite-correlation-matrices
X = mvnrnd(zeros(1,p), sigma, n);

% A is simple
A = normrnd(0,1,c, k);
% Simulate B as lower-trapezoidal with positive diagonal
B = normrnd(0,1,p,k); %actually full, we'll apply constraints later
% Apply shrinkage
B = B.*B_sig*G_sig;

% massage B here, so it will generate a problem solvable within constraints
% This is cheating, but I think I can use ML methods to cheat in the same
% way in the real problem.
ridx = zeros(1,k);
for col = 1:k
% Rearrange B so column col has the largest absolute value not already
% claimed by a previous column
tmp = B(:,col);
tmp(ridx(ridx~=0)) = 0;
[~, idx] = max(abs(tmp));
ridx(col) = idx;
end
new_idx = [ridx setdiff(1:p, ridx)];

%% re-arrange B and B_Sig to remain consistent
B = B(new_idx,:);
B_sig = B_sig(new_idx,:);

%% Apply constraints to B
for icol = 1:k
B(icol,icol) = abs(B(icol,icol));
B(1:(icol-1), icol) = 0;
end

noise = normrnd(0, noise_sig, n, c);
y = X*B*A'+noise;

%% Initial values seem to matter an awful lot.
%
%  The idea is to initialize the chains at some approximate maximum
%  likelihood estimate.

% get a maximum likelihood estimate of beta
beta = X\y;

% decompose int B*A';
[Q,L] = qr(beta');

% This produces iA with l-2 norm of 1. Actual A does not have norm = 1.
% When A is N(0,1), then expected norm is sqrt(numel(A)).
iB = L(1:k, :)'./sqrt(numel(A));
iA = Q(:,1:k).*sqrt(numel(A));

% flip iB so diags are positive
sgn = sign(iB(diag_idx))';
iB = iB.*sgn; %this will only work in recent versions of MATLAB.
iA = iA.*sgn;

% Add some jitter so each chain starts in a slightly different spot.
for i = 1:NC
idiags = abs(iB(diag_idx) + normrnd(0, init_sd, k,1));
ilowtrap = abs(iB(trap_idx) + normrnd(0, init_sd, sum(trap_idx(:)), 1));
ml_init(i) = struct('diags', idiags, 'lowtrap', ilowtrap, 'A', iA); %#okAGROW
end

Thanks again for everyone’s help so far, and I’m happy to get further feedback.

#10

For reference, the prior I mentioned in the context of the blog post is specifically designed to get around the column-ordering problems. It shows what prior a Gaussian prior is transformed to after applying a QR transform (i.e. transforming it to the trapezoidal form).

#11

Thanks for the reminder. I didn’t mean to blow past this part of the suggestion, but I’m having trouble figuring out how the logic in Leung & Drton maps onto this problem. I’ll try to briefly summarize my understanding of what they’re doing, in hopes it will expose some error in my thinking.

They are working with the form y=\beta f+e [1.1], where f is a k-vector of latent factors and \beta is m-by-k. The first thing they do is integrate out f, so now y\sim N(0,\Sigma), where \Sigma=\Omega+\beta \beta' [1.2]. This is great, because if you fix the order of y, you’ve uniquely specified \beta \beta', and if you shuffle the order of y, you can get the same answer by shuffling the rows and columns of \Sigma in the same way. But if we work backwards, re-ordering f has no effect on \beta \beta'. In latent factor analysis, that’s fine, because we want to be invariant to rearrangement of f. Here, though, we can think of the rows of XB as analogous to f, but I don’t think I want to integrate out XB to arrive at some covariance matrix with the order invariance I’m after.

The prior Leung & Drton explicate for the diagonals of L really makes the distribution of L into the LQ decomposition of \beta. What I’m doing is not that, and I probably should if I want to be fair to the model in simulations, but I don’t see how that would remove the need to be careful about the column order of X. I would love to be wrong about that, though.

Anyway, the prior they indicate is L_{i,i}\sim x^{k-i}exp\{-\frac{1}{2C_0}x^2\}. Since I’m modeling each entry in the lower trapezoidal part of B as coming from a normal distribution with (non-identical) variance, I don’t have a convenient C_0 for the distribution of the diagonal. I could back off from the assumption that I need diag_sig, lowtrap_sig and global_shrink in my model, or I could dig into the math and figure out how to adapt that prior. The former isn’t hard, so if it works, it might indicate that the latter option is worth pursuing.

#12

Let’s see if I am on the right track myself. The original non-identifiability arises as XB=(XO)(O^{-1}B) for any invertible O. Getting rid of the non-identifiability corresponds to finding a way to pick out O uniquely. The triangular constraint corresponds to picking the orthogonal matrix O=Q^{T} matching the LQ factorization of X, such that the equation becomes XB=L(QB). Now, if B has a uniform Gaussian distribution then \hat{B}=QB has the same distribution, and we can write XB=L\hat{B} with \hat{B} having the same prior as the original B. The problem now is to pick a prior on L . Any choice of prior then maps back to X=LQ^{T}, so many priors will induce a weird skewed distribution on the original X. Leung and Drton show that if X has a standard Gaussian prior, then L must have the particular prior they detail.

There is an issue, as you note, with imposing irregular variances: if B is not row-wise uniform (column-wise should be fine), then QB will not be identically distributed. Likewise, the Leung and Drton result only holds if X is i.i.d. in its elements. We can effectively change the row-wise variance of X, though, and column-wise variance of B by multiplying with scaling matrices: D_1 XB D_2. So we have a bit of a control still.

#13

Thanks again for taking the time to help me out. Earlier, I was thinking about B and A^T as analogues to the LQ factors of \beta, where Y=X\beta. In that context, I understand why I should use the Leung & Drton priors on the diagonals of L.

But now we’re talking about LQ factorization of the design/predictor matrix X. I can see how that would eliminate worries about the order of columns of X. So now that L and Q are specified by X, I don’t really know what it means for them to have priors. Bringing it back to Stan, I could do the LQ factorization of X in the transformed data block, but I’m not sure what it then means to put a sampling statement for L in the model block.

#14

Ah, sorry, I should have done my due diligence and double-checked the notation you used. In my notation X corresponds to your B and B to your A', i.e. I am solely considering the factor model part of the problem. For your convenience, here is my previous text with the necessary swaps:

edit: with respect to your comment, the key is not to match L and Q to B and A^{T}, but rather to decompose B using LQ decomposition.

#15

Thanks for the updated notation. The idea that B=LQ was indeed what I was missing. So at the end of the day, I don’t have a posterior sample on B and A^T, I have a sample from L and QA^T. I don’t see a direct way to recover B and A. Since that’s a matter of interpretation rather than modeling, I might just be satisfied with L and QA^T, or I could do some other decomposition after-the-fact.

Simulations with standard normal priors on B and A are working great, even with relatively wide X predictor matrices. To get a stable solution, it’s still important that the k diagonals of L are not too close to zero, and this can be achieved by permuting the columns of X. In practice, it looks like I can get a decent ordering by doing a QRP factorization on \beta^T, where \beta is a least-squares estimate in the form Y=X\beta. The ideal ordering would come from QRP factorization of B, but of course B is not easy to estimate.

I’m going to try to tackle the scaling matrix idea next, and assuming that works I’ll post the full model.

#16

Yeah, it is my impression that it will work with Gaussian priors as well (I assume you mean with the triangular constraints?), the issue is that it is quite hard to interpret the prior because of the ordering effect.

Since B is not really uniquely defined, what are you hoping to retrieve with your post hoc analysis? Which of the latents do you need for your further analysis?

#17

Right, I phrased that confusingly. In my simulations I select B and A from standard normal to generate the simulated data, so I model L as lower trapezoidal with the Leung & Drton priors on L.

I’ve actually realized that my proposed column re-ordering approach yields a distribution of the diagonals of L that is not what L&D describe. It makes sense; I’m manipulating the situation to try to get bigger values on the diagonals, so of course it will change the distribution. Maybe that’s what you meant about the prior being hard to interpret. Now, the L&D prior is still pretty reasonable, so maybe that’s not a problem.

Well, the first thing I thought to do was SVD on LQA^T, since that gets me something shaped like the output of partial least squares regression. I don’t really have a particular justification I can articulate for wanting B and A^T instead of L and QA^T. I think it’s just that the zeros in the upper triangle are harder for me to intuitively grasp. Nothing is usually exactly zero, but there those zeros are providing the uniqueness of the solution.

#18

LQA^T is good because it is an identifiable quantity. Doing SVD on each sample should then give you the distribution over singular vectors and singular values, but when I tried it I found it to not be as easy as I had hoped. First off, the SVD is also non-identifiable since you can do sign switches of each left and right vector simultaneously. This is a minor type of ambiguity, but quickly gets annoying to resolve once you take the uncertainty into account. Next, think about what happens when the eigenvalues are identical, then you can also apply simultaneous permutations to the vectors. While it is unlikely that you have perfectly equal eigenvalues, once you start looking across samples the ordering will become inconsistent as singular values switch spaces.

I agree completely with your misgivings about the forced triangular structure, but recall that due to the non-identifiability it is mostly a cosmetic issue. I feel that the nicest way to interpret it is by noting that if each column a of QA^T is Gaussian with \mathcal{N}(\mu,\Sigma), then x=La+e is Gaussian with x\sim\mathcal{N}(L\mu,L\Sigma L^T+D) for e\sim\mathcal{N}(0,D) and both a and e marginalized. So L\Sigma L^T is the low-rank covariance structure of each column of x, and it has no pesky zeroes or non-identifiabilities. The triangular bit of L is functionally just a cholesky decomposition of the covariance for the first few elements.

#19

The model based on the logic laid out in this thread seems to work okay with some cursory checks on parameter recovery. However, when I do simulation-based calibration (SBC), it looks like there is some problem.

My understanding of SBC is that in a well-specified model, when the parameters are drawn from their respective prior, used to simulate some data and then the model is fit to that data, then the true parameter value should be somewhere in the posterior distribution with high probability of falling in regions with high mass and low probability of falling in regions with low mass. In other words, the rank of the parameter with respect to the posterior should follow a uniform over many such simulations.

When I do this for the model (code below), some of the parameters look fine, but others have rank distributions with spikes at zero and max, like this:

(the black line shows the expected value under uniformity, with the shaded region showing a 99% ci. The right panel shows a binomial CDF in red vs. the observed CDF)

The reference paper describes this pattern as indicative of high autocorrelation in the chains. I have attempted what seems like pretty extreme thinning (as high as thinning at 2048), but I can’t seem to get the number of effective samples up to the target (127), usually with one or more parameters ending up with estimated number of effective samples in the 100-115 range.

At this point, I’m assuming there is some problem with the model, or that I need to do some kind of reparameterization. I suppose it is possible that I need to increase the number of warmup draws per chain, but I’m already going up to 16k per chain (which means SBC takes a long time).

Here is the model:

functions {
real ld_diag_lpdf(real x, int i, int k, real c0) {
return (k-i)*log(x) - square(x)/(2*c0);
}
}
data {
int n;
int p;
int c;
int k;
matrix[n, p] X;
matrix[n, c] Y;
}
transformed data {
int ntrap = p*k - k*(k-1)/2 - k;
}
parameters {
real<lower=0> sigma;
real<lower=0> nu;
vector<lower=0>[k] diags;
vector[ntrap] lowtrap;
matrix[k, c] Ahat;
vector<lower=0>[p] lambda;
}
transformed parameters {
matrix[p,k] L;
{ //hide integer
int idx;
idx=0;
L = rep_matrix(0, p, k);
for (col in 1:k) {
L[col, col] = diags[col]; //already applied sig_B
for (r in (col+1):p) {
idx+=1;
L[r, col] = lowtrap[idx];
}
}
} //close hider
}
model {
matrix[n, c] mu;
lowtrap ~ normal(0,1);
for (i in 1:k)
diags[i] ~ ld_diag(i, k, 1);
to_vector(Ahat) ~ normal(0,1);
lambda ~ normal(0, 2);
nu ~ gamma(2, 0.1);
sigma ~ cauchy(0,1);
mu = diag_post_multiply(X, lambda)*L*Ahat;
for (i in 1:c)
Y[:,i] ~ student_t(nu, mu[:,i], sigma);
}

And here is the simulation code:

data {
int n;
int p;
int c;
int k;
}
model {
}
generated quantities {
// parameters
real<lower=0> sigma;
real<lower=0> nu;
vector<lower=0>[p] lambda;
matrix[p,k] B;
matrix[c,k] A;

// transformed parameters
matrix[k, c] Ahat;
matrix[p,k] L;
matrix[p,p] B_fat;
matrix[p,p] qt_fat;
matrix[p,p] lt_fat;

// model values
matrix[n, p] X;
matrix[n,c] Y;
matrix[n,c] mu;

// parameter samples
sigma = fabs(cauchy_rng(0,1));
nu = gamma_rng(2, 0.1);
for (i in 1:p)
lambda[i] = fabs(normal_rng(0, 2));
for (j in 1:k){
for (i in 1:p)
B[i,j] = normal_rng(0,1);
for (i in 1:c)
A[i,j] = normal_rng(0,1);
}

// transformations
B_fat = append_col(B, rep_matrix(rep_vector(0,p), p-k));
qt_fat = qr_Q(B_fat');
lt_fat = qr_R(B_fat');
L = lt_fat[1:k,]';
Ahat = qt_fat[1:k,1:k]' * A';

// model code
for (i in 1:n) {
for (j in 1:p) {
X[i,j] = normal_rng(0,1); // candidate for change, different from matlab version
}
}
mu = diag_post_multiply(X, lambda)*L*Ahat;
for (i in 1:n) {
for (j in 1:c) {
Y[i,j] = student_t_rng(nu, mu[i,j], sigma);
}
}
}

As an aside, I’m using the Leung & Drton prior for the diagonals of L via a custom lpdf function. I couldn’t figure out how to sample from that distribution, so I had to back off one step and simulate B and A and then derive L. This seems like the most likely reason for a mismatch between generation and model, but I haven’t managed to find an error there yet. Would it be more direct to use the built-in rectangular Cholesky factor functions somehow?

#20

Yes, that’s the right understanding of SBC and you’re right this looks like it’s failing. That could be either because of the mismatch you alluded to (I haven’t heard of those distributions)—you really do need to have exact matches in the simulation and the fitting. But it could also just be a problem with fitting that could be solved with reparameterization.

But it also looks like the model’s not converging, and this is all based on getting models to converge so they’re sampling from the posterior.

You can speed up that last line with

to_vector(Y) ~ student_t(nu, to_vector(mu), sigma);

The Cauchy priors can be problematic if there’s not enough data to control the posterior.


text:
# Practice on Toph

Participate in exhilarating programming contests, solve unique algorithm and data structure challenges and be a part of an awesome community.

# Love Letter

By Nafis_Shhriar · Limits 1s, 512 MB

Mr. A has finally found the love of his life. Being a love expert Mr. P has given Mr. A an idea to write a love letter. But he is very weak in this matter. So he took help from his friends. He has $2 \times N$ friends. Each friend has given him one romantic sentence and using all those sentences Mr. A will write a love letter.

Each sentence contains only lower case alphabet (‘a’- ‘z’) and space. The height of each character is 1 and the width of each character is also 1.

But Mr. A is very busy person, So Mr. A has given me the duty to arrange the sentences and make a love letter. He has also given me an $N \times M$ size paper. Now I have to write all the $2 \times N$ sentences in that paper. But Mr. P has advised me not to write more than two sentences in a row and write in only one side of the paper. He also advised me not to write any sentence using more than one row. So, in order to fit the sentences, I will remove all the spaces and also I will not use any space between two lines.

Now I'm in a dangerous situation, if I refuse to do this work Mr. A will break my legs but if I do the work and something goes wrong he will break my legs too. So I've decided to do this work and hide for a few days. As I'm busy packing my bag i want you to do a job for me.

You will be given N and M, the height and width of the paper, and $2 \times N$ sentences. You will have to tell me can I write all the sentences in that paper in any order or not.

## Input

The first line of the input set will contain two integers $N$ and $M$ ($1 \le N, M \le 1000$).

The next $2 \times N$ lines will contain a sentence each, not containing more than 1000 characters.

## Output

If I can write all these sentences in the paper in any order, output $\texttt{YES}$ and if I can’t output $\texttt{NO}$.

## Samples

InputOutput
2 6
a a a
b b b
c c c
d d d

YES

InputOutput
2 28
akasata meg ura
bondu amar balobasay

NO

### Statistics

88% Solution Ratio

abcdef123Earliest, Jul '20

abcdef123Fastest, 0.0s

format: plain


text:
Which forest biometrics researchers get the most citations?

November 30, 2019
faculty publishing Google Scholar forest measurements forest inventory citations web scraping productivity

Researchers in the forest and natural resource sciences are measured by their productivity: how many publications they write and how many times those publications are cited by other researchers. Researchers that have many publications that are cited in the literature by others generally have rewarding careers and are promoted.

In 2017 I wrote a blog post for the Ecological Society of America’s Early Career Ecologist section that provided some evidence that early-career scientists publish earlier in their careers compared to more senior researchers.

When I wrote that post I used data from my own Google Scholar profile. I manually went through the data, including adding together citation counts and looking up researcher profiles to see when their PhD’s were awarded. Needless to say, it turned out to be very messy for what I thought would be a quick analysis and blog post.

I was happy to see that a package existed that can scrape data from Google Scholar and provide a number of useful metrics.

R’s scholar package

Many faculty at colleges and universities have a Google Scholar page. (For an example, see mine.) Having a Google Scholar page allows faculty to display their research publications, track the number of citations they receive for their articles, and interact with their coauthors.

The scholar package was developed by James Keirstead in R to extract citation data from Google Scholar. You can also obtain information about a specific researcher, compare multiple researchers, and predict future h-index values. (An h-index is the largest number h in which h publications have a least h citations.)

You can install and call the scholar library by using:

#install.packages("scholar")
library(scholar)

The scholar package can be used by identifying the Google Scholar ID for each researcher. The Google Scholar ID can be found at the end of the URL on each profile page:

You can specify the Google Scholar ID for a researcher and then obtain information from their profile using the get_profile() function:

# Define the ID for Matt Russell
id <- '4OO83OcAAAAJ'

# Get Matt's profile information
get_profile(id)
## $id ## [1] "4OO83OcAAAAJ" ## ##$name
## [1] "Matthew Russell"
##
## $affiliation ## [1] "University of Minnesota" ## ##$total_cites
## [1] 739
##
## $h_index ## [1] 16 ## ##$i10_index
## [1] 26
##
## $fields ## [1] "verified email at umn.edu - homepage" ## ##$homepage
## [1] "http://health.forestry.umn.edu/"
##
## \$coauthors
##  [1] "Christopher Woodall"  "Grant M. Domke"       "Anthony D'Amato"
##  [4] "Aaron Weiskittel"     "Shawn Fraver"         "Ram K. Deo"
##  [7] "John W Coulston"      "Charles Hobart Perry" "Laura Kenefic"
## [10] "Warren Cohen"         "Hans-Erik Andersen"   "Joshua Puhlick"
## [13] "Harold Burkhart"      "Michael J. Falkowski" "John A. Kershaw, Jr."
## [16] "Brian Palik"          "Tuomas Aakala"        "Mark Ducey"
## [19] "Ken Skog"             "sassan saatchi"

Then, you can obtain the number of citations by year for any researcher. To do this, you’ll use the get_citation_history() function:

# Get Matt's citation history, by year
get_citation_history(id)
##    year cites
## 1  2010     3
## 2  2011     3
## 3  2012    17
## 4  2013    18
## 5  2014    31
## 6  2015    62
## 7  2016   113
## 8  2017   141
## 9  2018   170
## 10 2019   164

Forest biometrics faculty productivity

I was interested in looking at the citations for some authors I was familiar with. Specifically, I investigated all authors from the following textbooks in the disciplines of forest measurements, growth and yield, and stand dynamics. These textbooks are commonly used in undergraduate and graduate courses across forest and natural resource departments:

• Forest Measurements, Sixth Ed. by Burkhart et al.
• Forest Growth and Yield Modeling by Weiskittel et al.
• Forest Dynamics, Growth and Yield by Pretzsch
• Modeling Forest Trees and Stands by Burkhart and Tome
• Forest Mensuration, Fifth Ed. by Kershaw et al.
• Forest Stand Dynamics by Oliver and Larson
• Forest Management: To Sustain Ecological, Economic, and Social Values by Davis et al.

Including all coauthors on these books amounts to 18 unique researchers. After searching for each author, 10 of them maintain a Google Scholar profile. The following data set includes their name, Google Scholar ID, number of citations throughout their career, and the number years in their career.

## # A tibble: 10 x 4
## # Groups:   name [10]
##    name                    id           num_citations career_years
##    <fct>                   <fct>                <int>        <int>
##  1 Hans Pretzsch           _HQsysQAAAAJ         13863           21
##  2 Harold Burkhart         mqzZk3gAAAAJ         10945           39
##  3 Jerome Vanclay          N1MKi-oAAAAJ         10536           28
##  4 Bruce Larson            ookqs8gAAAAJ         10102           33
##  5 Chadwick Dearing Oliver MkQKfx0AAAAJ          9285           34
##  6 Margarida Tome          GR9rjA4AAAAJ          5540           23
##  7 Pete Bettinger          rvCOS60AAAAJ          5430           30
##  8 John A. Kershaw, Jr.    Ah6KRK4AAAAJ          4235           39
##  9 Mark Ducey              5JseyAcAAAAJ          3290           19
## 10 Aaron Weiskittel        e7z1aLEAAAAJ          2891           14

For these authors, the correlation between the length of their career and total number of citations is 0.31:

Visualizing biometrics faculty productivity

We can extract all of the Google Scholar IDs for the ten researchers by using a concatenate function c().Then we can obtain a data frame comparing the number of citations of each researcher to their work in a specific year. This is accomplished with the compare_scholars() function.

Finally, we can compare career trajectories by using the compare_scholar_careers() function:

# Define the IDs for all researchers
ids <- c('e7z1aLEAAAAJ','mqzZk3gAAAAJ','Ah6KRK4AAAAJ','N1MKi-oAAAAJ',
'_HQsysQAAAAJ','GR9rjA4AAAAJ','5JseyAcAAAAJ','MkQKfx0AAAAJ',
'ookqs8gAAAAJ','rvCOS60AAAAJ')

# Compare number of citations in a specific year
compare_scholars(ids)

# Compare researchers career trajectories, based on year of first citation
citations<-compare_scholar_careers(ids)

We can compare researcher citation counts across a number of years. Notably, early- and mid-career faculty have a shorter number of years compared to more senior researchers. This is also generally reflected in the total number of citations:

We can also “standardize” the number of citations for each researcher because Google Scholar provides these values for each researcher. The following graph shows the number of citations at years 5, 10, 20, and 30. Note that some early- and mid-career researchers have not been in their careers for more than 20 years, so no data exist for them at those time points:

Predicting future faculty productivity

Another function in the scholar package is predict_h_index(), a function that allows you to predict a researcher’s h-index up to ten years ahead into the future. The algorithm is based in part on an article published in Nature in 2012:

# Predict Matt's future h-index
ids <- '4OO83OcAAAAJ'
predict_h_index(ids)
##    years_ahead  h_index
## 1            0 16.00000
## 2            1 19.74206
## 3            2 22.98970
## 4            3 26.01435
## 5            4 28.72624
## 6            5 31.39877
## 7            6 34.43822
## 8            7 37.27784
## 9            8 40.61042
## 10           9 44.45729
## 11          10 47.64782

We can graph the future h-index values for our forest biometrics researchers for ten years ahead into the future. The predict_h_index() function shows the change in h-index values through 2029:

I’m not sure that forest scientists will commonly receive h-index values greater than 100 throughout their careers. So I have some reservations about its application in the forest sciences. It is interesting to note the non-linear increase in the projected h-index values.

Limitations of Google Scholar and scholar package

There are a few limitations of using an approach like this, and in using the scholar package. Most notably:

• Not everyone has a Google Scholar page. If you list all coauthors for the seven books mentioned above, 44% of authors do not have a Google Scholar page. I know that at least two of them are deceased…
• Predicting future faculty productivity is based on neuroscience. The h-index forecasts published in Nature were developed for the field of neuroscience. With the rising use and ease of access to citation data from sites like Google Scholar, similar research productivity forecasting tools could be developed that are tailored to the forest science community.
• There are many more authors in forest biometrics. This was just a sample of a few researchers, and many more could be used for a similar analysis.

Conclusion

Forest scientists and researchers (and the supervisors and department heads that evaluate them for promotion) commonly use Google Scholar to measure their research productivity. The scholar package in R is an excellent tool for assessing an individual’s research productivity and comparing research output to others. Metrics can be queried from Google Scholar and appropriate visualizations produced to aid in interpreting citation data.

By Matt Russell. Leave a comment below or email Matt with any questions or comments.

How to calculate how much carbon is in a tree

October 3, 2019
forest measurements regression carbon

A 50-year history of forest measurements by analyzing Forest Measurements

September 15, 2019
forest inventory sampling books forest measurements

Quickly explore your timber cruise data with six lines of code

The inspectdf package in R allows you to quickly summarize timber cruise data by species, tree status, and other variables.
Data viz timber cruise forest inventory inspectdf Communicating data


text:
Constructive mathematics

Idea

Broadly speaking, constructive mathematics is mathematics done without the principle of excluded middle (or other principles, such as the full axiom of choice, that imply it). Sometimes one adds further restrictions or alternatively adds axioms that contradict excluded middle but are otherwise consistent; these variations may be seen in the list of schools below. However, it is probably best to use more precise terms (predicativism, intuitionism, etc) in this case.

Constructivism is the philosophy that such mathematics is useful, or (more strongly) that non-constructive mathematics is wrong. Historically, constructive mathematics was first pursued explicitly by mathematicians who believed the latter. However, many modern mathematicians who do constructive mathematics do it not because of any philosophical belief about the wrongness of non-constructive mathematics, but because constructive mathematics is interesting in its own right. In the ‘pluralist’ approach to the foundations of mathematics, a constructive proof (when it exists) is better because it is valid in more versions of mathematics, but a classical proof remains valid for classical mathematics.

Another motivation for modern mathematicians —especially category theorists like those on the nLab— is that the study of constructive mathematics has potential applications to non-constructive mathematics. For example, even if one believes the principle of excluded middle to be true, the “internal” version of excluded middle in many interesting categories is still false; thus constructive mathematics can be useful in the study of such categories, even if mathematics is “globally” non-constructive. This is the neutral motivation for constructive mathematics from the nPOV.

Here we write mostly about the mathematics, trying to be mostly neutral philosophically.

Origins and schools

During the “foundational crisis” in mathematics around the beginning of the 20th century, a number of mathematicians espoused philosophies that are generally grouped together and labeled constructivist. The common feature of these philosophies was a rejection of axioms and principles of logic that lead to nonconstructive proofs. There was much talk at the time about potential vs absolute infinity, but from an axiomatic perspective, it turns out that (if one stops short of finitism) the two main culprits are

Thus, constructivists (including many still active today) reject proofs that make use of either of these. (In fact, it was realized in 1975 by Diaconescu that the axiom of choice implies the principle of excluded middle. Excluded middle is precisely the “finitely indexed” axiom of choice; see excluded middle for a proof.)

There are, however, differences among constructivists as well.

• Some, like Errett Bishop, simply remove choice and excluded middle from classical mathematics with nothing to replace them. However, this makes it difficult to define a satisfactory notion of continuous function, even $ℝ\to ℝ$ without using locales; see (Waaldijk 2003).

• Others add “non-classical” axioms which contradict choice or excluded middle, but which are consistent in their absence, such as “all total functions $\left[0,1\right]\to ℝ$ are continuous” (the continuity axiom? of the “intuitionistic” school of L. E. J. Brouwer) or “all partial functions $ℕ\to ℕ$ are computable” (the computability axiom? of the “Russian” school of A. A. Markov, which is also called “constructive recursive analysis”).

• Most accept weaker versions of choice, such as countable choice or even dependent choice, but the school of Fred Richman rejects even these; see (Richman 2000). Toby Bartels argues that the intuition behind accepting these justifies the (stronger) axiom COSHEP (the category of sets has enough projectives).

• Still others, following Hermann Weyl, go even further and refuse to allow “impredicative” constructions; see predicative mathematics. Most of the work done by the schools of Brouwer, Bishop, and Markov (but not Richman) is also predicative, even though those founders were not adamant about it; as a result, predicativism often appears in the foundations of constructive mathematics (in particular, in those of Aczel and Martin-Löf, but not those of Friedman or Coquand).

• Very extreme constructivists, like Doron Zeilberger, reject the existence of infinite sets; see finitism. (Kronecker, the ‘father of constructivism’, is often considered a finitist, but as he came before the foundational crisis, it is difficult to classify him accurately.) Ironically, this means that excluded middle and choice, in their naïve forms, may become acceptable, or even true. (Even constructivists believe that they are true in the topos FinSet of finite sets.)

• The most extreme position of all is to deny even the existence of very large finite sets. This is called ultrafinitism. According to Wikipedia (to be fair, in an edit originally made by our own Toby Bartels), “even constructivists generally view [ultrafinitism] as unworkably extreme.”

Many constructivists (like many classical mathematicians) believe in an absolute mathematical sense of “truth,” and that in this sense choice and excluded middle are simply wrong. (Some constructivists, using classically false axioms, can even refute them; others merely claim that no possible correct reasoning could ever prove them. See Truth versus assertability below.) To most mathematicians, this makes them seem quite strange. Other constructivists adopt a wait-and-see attitude, or even a relative notion of truth (which can seem strange in another way).

Topos theory

With the invention of topos theory in the second half of the 20th century, a new sort of constructivism arose. It was observed (by Lawvere and others) that any topos with a natural numbers object has an internal logic which is powerful enough to interpret most of mathematics, but that this logic in general fails to satisfy choice and excluded middle. This meant that even for a mathematican who likes to use choice and excluded middle (and a fortiori for one who believes them to be “true”), there is a reason to care about what can be proven without them, because only if a proof is constructive can it be interpreted in an arbitrary topos with NNO.

Even starting from a “completely classical” world of mathematics, many interesting toposes arise naturally (such as the category of sheaves on any topological space, or more generally any site) whose internal logic is not classical logic. Then by internal reasoning in such a topos (which is constructive reasoning), one can prove various useful facts, which can then be reinterpreted as external statements about the behavior of the topos itself. (Some examples would be nice!)

By now it is known that many of the non-classical axioms used by the early constructivists have natural models in particular toposes. For instance, in the topos of sheaves on the real numbers $R$, it is true that “all total functions $R\to R$ are continuous.” And in the effective topos, it is true that “all partial functions $N\to N$ are computable.”

However, there are no non-classical (or classical) axioms beyond “pure constructivism” that are true in all toposes with NNO. In particular, there is a free topos? with NNO such that a statement is true in the free topos precisely when it is provable in pure (Richman-school) constructive mathematics. This means that for an argument to apply in all toposes, even mild assumptions such as countable or dependent choice are unacceptable. However, topos theory has also provided ideas that solve many of the problems with pure constructivism. For example, a well-behaved notion of “continuous function” can be recovered by using locales instead of topological spaces.

Some features of constructive mathematics

Rephrasing of classical ideas

Sometimes, all that is necessary to make a piece of mathematics constructive is careful use of language. It is common in classical mathematics to define things with an unnecessary amount of negation. This doesn’t work so well constructively, since a statement can be not false without being true, but we can sometimes do perfectly well by just removing unnecessary double negations.

For example, classically one often speaks about “nonempty” sets, meaning a set which “does not have no elements.” Constructively it is much better to say that a set “does have at least one element”; constructivists often call such a set inhabited or occupied to remind themselves that it is a “positive” notion to replace the negative one of “nonempty”. Others continue to use the word “non-empty” but understand it as a term of art that really means “inhabited”.

Bifurcation of notions

On the other hand, differences in axiomatization or definition that make no difference classically can result in actual differences in behavior constructively. Therefore, classically equivalent notions often “bifurcate” (or “trifurcate” or worse) into multiple inequivalent constructive ones. This tends to happen whenever a concept involves negation and, to a lesser degree, disjunction and existential quantification. In some cases there is a “correct” constructive version of the definition, although it may take some thought to uncover it, but in many cases more than one of the resulting concepts is important and useful. For example:

• There are multiple inequivalent constructive definitions of a field, because of the axioms “every nonzero element has an inverse” and $0\ne 1$.

• In pure constructive logic, the Dedekind real numbers and the Cauchy real numbers need no longer coincide. The Cauchy reals sit inside the Dedekind reals, but in general not every Dedekind real need be approximable by a sequence of rationals. From a topos-theoretic viewpoint, the Dedekind reals are usually the “correct” notion to study (if not the locale of real numbers as a whole). However, either excluded middle or countable choice suffices to ensure (via weak countable choice) that every Dedekind real is a Cauchy real, and hence that the two notions coincide; see (Bridges et al, 1998).

• There are at least three different constructive notions of ordinal number; see (Taylor 1996) and (Joyal–Moerdijk 1995).

• Without the axiom of choice, functors and anafunctors become distinct, and often the latter is more appropriate; see (Makkai, 1996).

• Perhaps most disturbingly of all to the classical mathematician, one must distinguish between finite sets, subfinite sets, finitely-indexed sets, and even subfinitely indexed sets; see finite set for definitions. However, in practice it is usually either finite or finitely-indexed sets that are important, and with practice a little bit of thought suffices to show which is the relevant concept.

Negative translation

This allows one to translate classically valid theorems into intuitionistically valid theorems. See double negation translation.

Truth versus assertability

Already in classical mathematics, there is a difference between saying that something is true and saying that something is provable. If you adopt ZFC because you believe it correct, then (assuming that you're aware of certain theorems) you believe that ZFC is consistent even though you also know that you cannot prove it so. In that case, you also believe that the continuum hypothesis (for example) is either true or false; you may or may not have an opinion on which it is, but in any case again you know that you cannot prove it either way.

A constructive mathematician can be even subtler. If you belong to the Bishop school, then you accept no classically false axioms, and anything that you can prove is valid also in classical mathematics. Even so, you can believe that the principle of excluded middle is false (even though, of course, you know that you cannot prove it false). So you can really confuse the classical mathematicians by saying, on the one hand, that they can safely accept all of your theorems as valid, and then saying, on the other hand, that you know excluded middle to be false. The resolution, of course, is that you never claimed that this was a theorem.

This way of talking has even been formalised in (Bishop, 1967) with a convention adopted by many (but not all) of his followers. In this convention, the word ‘not’, used normally in a vernacular sentence whose mathematical content would otherwise be the proposition $p$, changes to the content to $¬p$, or $p\to \perp$, as usual. (This follows all of the rules of intuitionistic logic; in particular, any statement that is ‘not’ true is false (just as in classical logic).) However, the word ‘not’, in italics as shown here, changes the content to $p\to x$, where $x$ is some statement that is known (although not proved!), according to Bishop, to be false. Now a statement that is ‘not’ true may be false, but this may be unknown, and it is even possible that it is also ‘not’ false (so $¬p\to x$, possibly for a different $x$).

Bishop gives in his introduction several statements, suitable for use as $x$ above, including:

• excluded middle itself, which Bishop call the ‘principle of omniscience’;
• the ‘limited principle of omniscience’: any infinite sequence in $\left\{0,1{\right\}}^{N}$ is either all $0$ or has at least one $1$;
• the ‘lesser limited principle of omniscience’: for any two infinite sequence in $\left\{0,1{\right\}}^{N}$ that do not both have at least one $1$, at least one of these sequences does not have at least one $1$;
• others.

At one point in his book, while discussing the possibility of a pointwise-continuous function $\left[0,1\right]\to R$ that is not uniformly continuous (a ‘Specker function’, whose existence Markov's school claims as a theorem), Bishop seems to assert that this theorem is both not true and not false; he does not put it this way, so this may not be exactly what he meant, but there is no contradiction if it is. (But it is a contradiction, even in intuitionistic logic, if a statment is both not true and not false; indeed, a definition of ‘false’ may be taken to be ‘not true’.)

This practice can be understood through a careful distinction between object language and metalanguage. A mathematical statement $p$ (such as the continuum hypothesis, or that a Specker function exists) belongs to the object language, as does $¬p$ or $p\to x$ for any specific mathematical statement $x$. But the statement that $p$ is provable in some formal system belongs to the metalanguage (although it can also phrased internal to any object language suitable for mathematics), and as such may be written $⊢p$. The metalanguage has its own logic (the metalogic, which for the sake of argument we may even assume to be classical), but notice that $¬⊢p$ and $⊢¬p$ are different; the first claims that $p$ cannot be proved, while the second claims that $p$ can be refuted, which (in any consistent? formal system, even a classical one) is stronger. Although Bishop did not commit to any formal system, if we assume for the sake of argument that he settles on one, then we may write $⊢¬p$ as our interpretation of his meaning when he asserts ‘not $p$’ but $¬⊢p$ as his meaning when he asserts ‘not $p$’.

Even without this convention, which is peculiar to the Bishop school, it is important when reading constructive mathematics to remember the difference between what can be refuted (proved false) and what merely cannot be proved true. Although Brouwer's and Markov's schools will sometimes claim to prove statements that are (classically) outright false (such as that every function $\left[0,1\right]\to R$ is pointwise-continuous), it is more common (and can happen in any school) that a constructive mathematician makes a claim that merely sounds false, when what they really mean is only that they don't accept what you say as true. Often modal phrases like ‘not necessarily’ will be used where Bishop would use ‘not’, as a clue that we're shifting to a metalanguage, or at least merely remaining agnostic rather than outright disagreeing.

The distinction between object language and metalanguage exists even in classical mathematics, but it seems that most classical mathematicians are not used to remembering it, although it is not entirely clear why this should be so. One possibly relevant observation is that even if $P$ is a statement which is neither provable nor disprovable (like the continuum hypothesis), in classical mathematics it is still provable that “$P$ is either true or false.” Moreover, classical model theory often restricts itself to two-valued model?s in which the only truth values are “true” and “false,” although classical logic still holds in Boolean-valued models, and in such a case the truth value of $P$ may be neither “true” nor “false,” although the truth value of ”$P$ or not $P$” is still “true.” Certainly when talking about classical truths which fail constructively, such as excluded middle, it is important to remember that “fails to be provable” is different from “is provably false.”

Concepts that usually arise in constructive mathematics, often because they are classically trivial:

Some of these are also useful internally or even classically.

Topics relevant to the foundations of constructive mathematics:

Other articles with content relating to constructive mathematics:

In principle, every article could explain how it applies to constructive mathematics, but that will probably never happen.

References

Most books on topos theory include some discussion of toposes' internal constructive logic. One good reference is:

This discussion originally took place at apartness relation:

Revised on April 23, 2013 20:21:56 by Urs Schreiber (131.174.42.75)


text:
# Not sure how to start this problem hint please!

An empty rubber balloon has a mass of 0.0135 kg. The balloon is filled with helium at a density of 0.181 kg/m3. At this density the balloon is spherical with a radius of 0.300 m. If the filled balloon is fastened to a vertical line, what is the tension in the line?

I don't even know what equation to use...! AHH

Related Introductory Physics Homework Help News on Phys.org
Tide
Homework Helper
HINT: It's a question about buoyancy.

is the tension equal to the force?

lightgrav
Homework Helper
yes, but there are several forces involved here. ...
all dealing with one Volume.

Tide
Homework Helper
Yes. The buoyant force is balanced by the tension in the string.

so would i use the equation B=density(fluid) x volume(fluid) x g?

Pyrrhus
Homework Helper
Yes, substitute it in your $\sum F_{y} = 0$

so it would be Density(Helium) x what volume? x 9.81 ?

mukundpa
Homework Helper
The forces acting on the balloon are gravity force, buoyant force and the tension in the line.

The buoyant force is applied by air of the atmosphere.

Buoyant force is equal to the weight of the displaced volume of the fluid.


text:
# Show that the area of the triangle contained between the vectors

Question.
Show that the area of the triangle contained between the vectors $\mathbf{a}$ and $\mathbf{b}$ is one half of the magnitude of $\mathbf{a} \times \mathbf{b}$.

solution:

Consider two vectors $\overrightarrow{\mathrm{OK}}=|\vec{a}|$ and $\overrightarrow{\mathrm{OM}}=|\vec{b}|$, inclined at an angle $\theta$, as shown in the following figure.

In $\Delta \mathrm{OMN}$, we can write the relation:

$\sin \theta=\frac{\mathrm{MN}}{\mathrm{OM}}=\frac{\mathrm{MN}}{|\vec{b}|}$

$\mathrm{MN}=|\vec{b}| \sin \theta$

$|\vec{a} \times \vec{a}|=|\vec{a}||\vec{b}| \sin \theta$

$=\mathrm{OK} \cdot \mathrm{MN} \times \frac{2}{2}$

$=2 \times$ Area of $\Delta$ OMK

$=2 \times$ Area of $\Delta \mathrm{OMK}$

$\therefore$ Area of $\triangle \mathrm{OMK}=\frac{1}{2}|\vec{a} \times \vec{b}|$


text:
# Get Answers to all your Questions

#### Write the correct answer :Which of the following is different from the others?(a) 20 + ( –25)  (b) (– 37) – (– 32)   (c) (– 5) × ( –1)   (d) ( 45 ) ÷ ( – 9)

(c)

$\\ (a) \ 20+(-25)=20-25=-5 \\ (b) \ (-37)-(-32)=-37+32=-5 \\ (c) \ (-5) \times(-1)=5 \\ (d) \ (45) \div(-9)=-\frac{45}{9}=-5 \\$

(-5)  x  (-1) is different from other options.


text:
#### Question

Consider a rigid plane body or "lamina," such as a flat piece of sheet metal, rotating about a point O in the body. If we choose axes so that the lamina lies in the xy plane, which elements of the inertia tensor I are automatically zero? Prove that

$I_{zz} = I_{xx} + I_{yy}.$

Verified

#### Step 1

1 of 2

We are considering a rigid plane body or “lamina,” such as a ﬂat piece of sheet metal that lies in xy plane.

Density of the metal sheet can be written as:

$\rho(x,y,z)=\sigma(x,y)\delta(z)$

If we integrate term above, we can get full mass of the sheet:

$\int{\sigma(x,y)\delta(z) dxdydz}=\int{\sigma(x,y)dxdy}\int{\delta(z)dz}=\int{\sigma(x,y)dxdy}=\mu$

Diagonal inertia tensors are given by:

$I_{xx}=\int{\sigma(x,y) (y^2+z^2)\delta(z)dxdydz}=\int{\sigma(x,y)\cdot y^2dxdy}$

$I_{yy}=\int{\sigma(x,y) (x^2+z^2)\delta(z)dxdydz}=\int{\sigma(x,y)\cdot x^2dxdy}$

$I_{zz}=\int{\sigma(x,y) (x^2+y^2)\delta(z)dxdydz}=\int{\sigma(x,y)\cdot (x^2+y^2)dxdy}$

We can see that :

$\boxed{I_{zz}=I_{xx}+I_{yy}}$

The only non diagonal element that is not equal to zero is:

$I_{xy}=-\int{\sigma(x,y) x\cdot y\delta(z)dxdydz}=-\int{\sigma(x,y)\cdot xy\cdot dxdy}$

All other non diagonal elements are equal to zero:

$\boxed{I_{xz}=I_{yz}=0}$

## Create an account to view solutions

#### Physics for Scientists and Engineers: A Strategic Approach with Modern Physics

4th EditionRandall D. Knight
3,509 solutions

#### Physics for Scientists and Engineers

9th EditionJohn W. Jewett, Raymond A. Serway
6,002 solutions

#### Physics for Scientists and Engineers with Modern Physics

10th EditionJohn W. Jewett, Raymond A. Serway
3,192 solutions

#### Classical Mechanics

1st EditionJohn R. Taylor
741 solutions


text:
warpPerspective Transform Matrix type ?

I am using opencv on C++. I am quite new at Opencv, and I think that my problem should be easy to solve ... thank you for your help. I don't understand what should be the data type of the transformation matrix in warpPerspective and how it works. I have a 3 x 3 matrix (of type Eigen::MatrixXf) H, and I manually build a cv::Mat H2 equal to H, so that I could execute that transformation on an image, as follows :

warpPerspective(src , newImg , H2 , Size (5*H1 , 5*W)  ) ;

As I have not succeeded to build H2 correctly, I do it and check it step by step, as follows :

Mat H2(3,3, CV_64F);
for (int g=0;g<3;g++) for (int g1=0; g1<3; g1++) H2.at<float>(g,g1)=  float(H(g,g1)) ;

cout <<endl<< "Matrix H2 (OpenCV) : "<<endl;
cout<< H2<<endl;
cout <<endl<< "Matrix H2 (OpenCV) Element by element: "<<endl;
for (int t=0; t<3; t++)
{
for (int s=0; s<3; s++)
{
cout<<"Element "<<t<<" , "<<s<<endl;
cout<< H2.at<Vec3f>(t,s)<<endl;
}
cout<<endl;
}
cout <<endl<< "Matrix H (Eigen) : "<<endl;
cout << H<<endl;

But when I print I get in H2 something that seems to be a 3x3x3 tensor, with column 0 for H2 equal to H. As shown in the attached screenshot. C:\fakepath\Capture du 2016-05-20 09:31:51.png

edit retag close merge delete

Sort by » oldest newest most voted

Mat H2(3,3, CV_64F); mean you have a matrix of double. Then you must use H2.at<double>(g,g1) (and not float)

If you want a matrix of float you must use Mat H2(3,3, CV_32F)

more


text:
# Using dominated convergence theorem to asset Leibniz rule when $\partial{f}/\partial{x}(x,y)$ is bounded and $f$ is integrable

I was thinking about the following problem

Let $f(x,y)$, $0 \leq x, y \leq 1$, be a function such that for each $x$, $f(x,y)$ is integrable, and $\partial{f}/\partial{x}(x,y)$ is bounded function of $(x,y)$. Show that $\partial{f}/\partial{x}(x,y)$ is a measurable function of $y$ for each $x$ and $$\frac{\partial}{\partial{x}}\int^{1}_{0} f(x,y)dy = \int^{1}_{0} \frac{\partial{f}}{\partial{x}}(x,y)dy$$ Assume that $\partial{f}/\partial{x}(x,y)$ exists everywhere for now (I'm not sure if we can assume this):

I think the problem is asking us to justify using Leibniz rule by applying dominated convergence theorem to the sequence of functions $$\phi_{h}(x,y) := \frac{f(x+h,y)-f(x,y)}{h}.$$ To do this, we need to show that the sequence is dominated by an integrable function, by the mean value theorem, there is an $x_{0} \in (x, x+h)$ where

$\phi_{h}(x,y)= \frac{\partial{f}}{\partial{x}}(x_{0},y)$

But I am a bit stuck relating $\frac{\partial{f}}{\partial{x}}(x_{0},y)$ to the integrable function $f$. I think I am missing some estimate.

You're overthinking this. Since $\partial{f}/\partial{x}(x,y)$ is bounded on all of $[0,1]\times[0,1]$, just let $M$ be a constant bound for it and apply the dominated convergence theorem with $M$ as your bound.

$\left\lvert\frac{\partial f}{\partial x}(\xi_h,y)\right\rvert\le \sup\limits_{\text{domain}}\left\lvert\frac{\partial f}{\partial x}\right\rvert=M$. Since constant functions are integrable on the finite interval $(0,1)$, by DCT you can assert that for all $x$ $$\lim_{h\to 0}\int_0^1\phi_h(x,y)\,dy=\int_0^1\lim_{h\to 0}\phi_h(x,y)\,dy$$


text:
## May 11, 2021 Single Round Match 804 Editorials

#### YetAnotherGraphColoring

Hint: Longest Decreasing Subsequence.

Fact 1: If for x, y x < y and a[x] > a[y] and for y, z y < z and a[y] > a[z] ⇒ x < z and a[x] > a[z].

Fact 2: Consider a sequence of indices i1, i2, i3, …, ik, where for each 1 < j ≤ k ij – 1 < ij and a[ij – 1] > a[ij]; then each pair of indices in the sequence are connected in our graph. So we need at least k colours.

What is the longest sequence satisfying the above conditions? Longest Decreasing Subsequence.

So if the length of LDS is k, we need at least k colours.

Read how to find the Longest Decreasing Subsequence here.

Just as when calculating LDS, let dpi be the length of longest decreasing sequence ending in i. For example for sequence 3, 2, 2, 1, 5, dp is [1, 2, 2, 3, 1].

If we assign colour dpi to i-th index, it’s proper colouring. Consider the opposite case, there exists x, y where x < y and a[x] > a[y] and dpx = dpy. As x < y and a[x] > a[y], we can add the y-th element to the LDS ending in x (with length dpx) and increase its length by 1, so dpy >= dpx + 1;  contradiction.

So the answer is the length of LDS. One can easily convert LDS to LIS by reversing the array.

Here is the code by at_f:

from collections import Counter, defaultdict, deque
from heapq import heapify, heappop, heappush
from bisect import bisect_left, bisect_right
import math
# import numpy as np

range = xrange

def lengthOfLIS(nums):
chain = [float("-inf")]
for x in nums:
newLen = bisect_left(chain, x)
assert chain[newLen - 1] < x
if newLen == len(chain):
chain.append(x)
else:
assert x <= chain[newLen]
chain[newLen] = x
return len(chain) - 1

class YetAnotherGraphColoring:
def minColors(self, n, a1, a2, x, y, z, m):
# Parameters: integer, integer, integer, integer, integer, integer, integer
# Returns: integer

a = [-1] * n
a[0] = a1
a[1] = a2
for i in range(2, n):
a[i] = (x * a[i - 1] + y * a[i - 2] + z) % m

lis = lengthOfLIS(a[::-1])
return lis

UnluckyNumber

Hint: Inclusion-exclusion.

We can count all of the possible options for assigning numbers, then subtract options where at least an edge is unlucky (i. e. sum of assigned numbers to its endpoints is z).

Let’s fix an edge to be unlucky, how many options are there for its endpoints? If we set the colour of one endpoint, the other has only one option. We can set one of the endpoints to max(1, z – k), 2, …, min(k, z – 1). So the number of options is min(k, z – 1) – max(1, z – k) + 1, let this number be len.

So we selected an edge and we want to calculate the number of colourings that this edge is unlucky in, the number of such colourings is k ^ (n – 2) * len. We have len options for endpoints of the corresponding edge and also k options for n – 2 remaining nodes.

This way, we subtracted some of invalid colourings twice. For example consider the case:

n = 4
m = 2
k = 1
z = 2
x[] = {1, 3}
y[] = {2, 4}

y[] = {2, 4}

The number of possible options must be 0, but with the above method we get -1.

We need to add the invalid colorings that counted twice and add them once. This results in a method known as inclusion-exclusion. Read about it here.

Let E be the set of all of the edges. Let f(A), where A is a subset of E, be the number of colorings where each edge in A is unlucky. The the answer is

Now the remaining part is to calculate f(A). The vertices don’t have any edge in A have simply k options. Others will form several components using edges in A. It’s easy to understand that if we determine the colour of a vertex in a component, other vertices in the component will uniquely get colour. What if we have cycles? There is no problem with even cycles. But odd cycles sometimes make problems. If we set colour x for vertex a which lies in an odd cycle, x must be equal to z – x, so x = z / 2 and z needs to be even in this case.

To summarize, each component without an odd cycle has len options. Components with odd length and odd z have 0 options. Components with odd length and even z have 1 option.

To calculate f(A) we multiply options for single vertices and components.

Time complexity is O(2^m*m).

Code by heno239:

using namespace std;

//#define int long long
typedef long long ll;

typedef unsigned long long ul;
typedef unsigned int ui;
constexpr ll mod = 1000000007;
const ll INF = mod * mod;
typedef pair<int, int>P;
#define stop char nyaa;cin>>nyaa;
#define rep(i,n) for(int i=0;i<n;i++)
#define per(i,n) for(int i=n-1;i>=0;i--)
#define Rep(i,sta,n) for(int i=sta;i<n;i++)
#define rep1(i,n) for(int i=1;i<=n;i++)
#define per1(i,n) for(int i=n;i>=1;i--)
#define Rep1(i,sta,n) for(int i=sta;i<=n;i++)
#define all(v) (v).begin(),(v).end()
typedef pair<ll, ll> LP;
typedef long double ld;
typedef pair<ld, ld> LDP;
const ld eps = 1e-12;
const ld pi = acosl(-1.0);

ll mod_pow(ll x, ll n, ll m = mod) {
if (n < 0) {
ll res = mod_pow(x, -n, m);
return mod_pow(res, m - 2, m);
}
if (abs(x) >= m)x %= m;
if (x < 0)x += m;
ll res = 1;
while (n) {
if (n & 1)res = res * x % m;
x = x * x % m; n >>= 1;
}
return res;
}
struct modint {
ll n;
modint() :n(0) { ; }
modint(ll m) :n(m) {
if (n >= mod)n %= mod;
else if (n < 0)n = (n % mod + mod) % mod;
}
operator int() { return n; }
};
bool operator==(modint a, modint b) { return a.n == b.n; }
modint operator+=(modint& a, modint b) { a.n += b.n; if (a.n >= mod)a.n -= mod; return a; }
modint operator-=(modint& a, modint b) { a.n -= b.n; if (a.n < 0)a.n += mod; return a; }
modint operator*=(modint& a, modint b) { a.n = ((ll)a.n * b.n) % mod; return a; }
modint operator+(modint a, modint b) { return a += b; }
modint operator-(modint a, modint b) { return a -= b; }
modint operator*(modint a, modint b) { return a *= b; }
modint operator^(modint a, ll n) {
if (n == 0)return modint(1);
modint res = (a * a) ^ (n / 2);
if (n % 2)res = res * a;
return res;
}

ll inv(ll a, ll p) {
return (a == 1 ? 1 : (1 - p * inv(p % a, a)) / a + p);
}
modint operator/(modint a, modint b) { return a * modint(inv(b, mod)); }
modint operator/=(modint& a, modint b) { a = a / b; return a; }
const int max_n = 1 << 1;
modint fact[max_n], factinv[max_n];
void init_f() {
fact[0] = modint(1);
for (int i = 0; i < max_n - 1; i++) {
fact[i + 1] = fact[i] * modint(i + 1);
}
factinv[max_n - 1] = modint(1) / fact[max_n - 1];
for (int i = max_n - 2; i >= 0; i--) {
factinv[i] = factinv[i + 1] * modint(i + 1);
}
}
modint comb(int a, int b) {
if (a < 0 || b < 0 || a < b)return 0;
return fact[a] * factinv[b] * factinv[a - b];
}
modint combP(int a, int b) {
if (a < 0 || b < 0 || a < b)return 0;
return fact[a] * factinv[a - b];
}

class UnluckyNumber {
public:
int numberOfColourings(int n, int m, int k, int z, vector <int> x, vector <int> y) {
vector<vector<int>> G(n);
rep(i, m) {
x[i]--; y[i]--;
}
int ri = min(z - 1, k);
int le = max(1, z - k);
modint len = ri - le + 1;
vector<bool> used(n);
vector<int> col(n);
modint ans = 0;
rep(i, (1 << m)) {
rep(j, n)G[j].clear();
fill(all(used), false);
rep(j, m) {
if (i & (1 << j)) {
G[x[j]].push_back(y[j]);
G[y[j]].push_back(x[j]);
}
}
modint num = 1;
rep(j, n) {
if (used[j])continue;
int cnt = 1; used[j] = true;
queue<int> q; q.push(j);
bool isb = true;
col[j] = 0;
while (!q.empty()) {
int id = q.front(); q.pop();
for (int to : G[id]){
if (!used[to]) {
cnt++;
used[to] = true;
col[to] = col[id] ^ 1;
q.push(to);
}
else {
if (col[to] != (col[id] ^ 1))isb = false;
}
}
}
if (cnt == 1)num *= k;
else {
if (isb) {
num *= len;
}
else {
if (z % 2 == 0)num *= 1;
else num *= 0;
}
}
}
int coef = 1;
rep(j, m)if (i & (1 << j))coef *= -1;
ans += (modint)coef * num;
}
return ans;
}
};

#### PolygonAndPermutations

Hint: Count inversions.

There are some intersections that it’s impossible to avoid. Consider the intersection between two segments (a, b), (b, d). If a, b, c, d belong to 3 or 4 different line segments, it’s impossible to avoid this intersection.

When a, b, c, d belong to 4 different line segments:

It’s not depend on the order we put the line segments or direction of segments (upside down or not), we’ll always have this intersection.

When a, b, c, d belong to 3 different line segments:

One may wonder if we make the top line segment upside down this intersection will resolve. But We’ll then have this intersection again:

So the only intersections we can minimize are those that a, b, c, d belong to two segments. So the order we put segments is not important and the only point is which segments we’re using “upside down”. For each two segments, we’ll calculate the number of intersections if we use them (1) both normally or both upside down and (2) one normally and one upside down.

Let’s calculate for two segments (with fixed direction) how many intersections are there. Calculation is easy in O(n^2); although one can calculate it in O(n log n).

At last we bruteforce on the direction of segments. Which segments we’ll use in upside down mode and which normal mode. Then we’ll calculate the number of intersections using the previously preprocessed values.

Code by tourist:

class PolygonAndPermutations {
public:
long long solve(int, int, long long);
};

long long xx, pr;

long long get_random() {
xx = xx % 1000000007;
long long temp = (xx*xx + pr) % 1000000007;
pr = xx;
xx = temp;
return xx;
}

vector<int> permutation(int n) {
vector<int> cur(n), res;
for(int i=0; i<n; ++i) cur[i] = i+1;
for(int i=0; i<n; ++i) {
int tmp = get_random() % (n-i);
swap( cur[n-i-1], cur[tmp] );
res.push_back(cur[n-i-1]);
}
return res;
}

long long PolygonAndPermutations::solve(int K, int N, long long seed) {
xx = seed;
pr = 0;
vector< vector<int> > segments(K);
for (int i=0; i<K; ++i) {
segments[i] = permutation(N);
assert((int) segments[i].size() == N);
for (int j = 0; j < N; j++) segments[i][j] -= 1;
}
for (int i = 0; i < K; i++) {
for (int j = i + 1; j < K; j++) {
for (int i2 = 0; i2 < K; i2++) {
for (int j2 = i2 + 1; j2 < K; j2++) {
if (make_pair(i, j) < make_pair(i2, j2)) {
if (i != i2 && i != j2 && j != i2 && j != j2) {
if (i < i2 && i2 < j && j < j2) {
}
continue;
}
assert(i == i2 || i == j2 || j == i2 || j == j2);
add += N * (N - 1) / 2;
}
}
}
}
}
vector<vector<int>> mat(K, vector<int>(K));
for (int i = 0; i < K; i++) {
for (int j = i + 1; j < K; j++) {
vector<int> pos_j(N);
for (int t = 0; t < N; t++) {
pos_j[segments[j][t]] = t;
}
for (int t0 = 0; t0 < N; t0++) {
for (int t1 = t0 + 1; t1 < N; t1++) {
if (pos_j[segments[i][t0]] < pos_j[segments[i][t1]]) {
mat[i][j] -= 1;
} else {
mat[i][j] += 1;
}
}
}
}
}
long long ans = (long long) 1e18;
for (int t = 0; t < (1 << K); t++) {
long long cur = 0;
for (int i = 0; i < K; i++) {
for (int j = i + 1; j < K; j++) {
if ((t >> i) % 2 != (t >> j) % 2) {
cur += mat[i][j];
}
}
}
ans = min(ans, cur);
}
}

#### CostMaximizer

Hint: I wish it was possible to select all vertices from the same color.

Consider we select a set of vertices S from blue and others from red. Then answer is

So we want to minimize

Spend O(n^2) and fix v, u such that v is a vertex with red color and u is a vertex from blue colour (we know that these vertices exists).

So we want to put vertices in two groups, one with v and one with u, such that the sum of weights of edges (2 * w – k) connecting vertices from different groups becomes minimum. This is the definition of minimum cut which is equal to maximum flow. Using max-flow we can solve the problem.

Code by fetetriste:

public class CostMaximizer {
final int infinity = (int) 1e9;

public long maxPossible(int n, int K, int[] inputW) {
int ptr = 0;
int[][] w = new int[n][n];
long sum = 0;
for (int i = 0; i < n; i++) {
for (int j = i + 1; j < n; j++) {
w[i][j] = inputW[ptr];
w[j][i] = inputW[ptr];
++ptr;
sum += w[i][j];
}
}

long min = infinity;
for (int s = 0; s < n; s++) {
for (int t = s + 1; t < n; t++) {
Graph g = new Graph(n + 2);
for (int i = 0; i < n; i++) {
for (int j = i + 1; j < n; j++) {
g.addEdge(i, j, 2 * w[i][j] - K);
g.addEdge(j, i, 2 * w[i][j] - K);
}
}
min = Math.min(min, g.minCut());
}
}

return sum - min;
}

class Graph {
int n;
int s;
int t;
int[] firstEdge;
int[] nextEdge;
int[] edgeDst;
int[] edgeCap;
int[] edgeFlow;
int numEdges;
int[] level;
int[] q;
int[] lastEdge;

Graph(int n) {
this.n = n;
s = n - 1;
t = n - 2;
firstEdge = new int[n];
final int initialEdges = 10;
nextEdge = new int[initialEdges];
edgeDst = new int[initialEdges];
edgeCap = new int[initialEdges];
edgeFlow = new int[initialEdges];
numEdges = 0;
level = new int[n];
q = new int[n];
lastEdge = new int[n];
Arrays.fill(firstEdge, -1);
}

void reallocate() {
int k = nextEdge.length * 3 / 2 + 1;
nextEdge = Arrays.copyOf(nextEdge, k);
edgeDst = Arrays.copyOf(edgeDst, k);
edgeCap = Arrays.copyOf(edgeCap, k);
edgeFlow = Arrays.copyOf(edgeFlow, k);
}

void addEdge(int a, int b, int cap) {
if (numEdges + 2 >= nextEdge.length) {
reallocate();
}
int e = numEdges++;
nextEdge[e] = firstEdge[a];
firstEdge[a] = e;
edgeDst[e] = b;
edgeCap[e] = cap;
edgeFlow[e] = 0;

e = numEdges++;
nextEdge[e] = firstEdge[b];
firstEdge[b] = e;
edgeDst[e] = a;
edgeCap[e] = 0;
edgeFlow[e] = 0;
}

long minCut() {
return maxFlow();
}

long maxFlow() {
long flow = 0;
while (bfs()) {
for (int i = 0; i < n; i++) {
lastEdge[i] = firstEdge[i];
}
while (true) {
break;
}
}
}
return flow;
}

int dfs(int v, int delta) {
if (v == t) {
return delta;
}
for (; lastEdge[v] >= 0; lastEdge[v] = nextEdge[lastEdge[v]]) {
int e = lastEdge[v];
int u = edgeDst[e];
if (level[u] != level[v] + 1) {
continue;
}
int x = Math.min(delta, edgeCap[e] - edgeFlow[e]);
if (x <= 0) {
continue;
}
x = dfs(u, x);
if (x > 0) {
edgeFlow[e] += x;
edgeFlow[e ^ 1] -= x;
return x;
}
}
level[v] = -1;
return 0;
}

boolean bfs() {
Arrays.fill(level, -1);
q[0] = s;
level[s] = 0;
int qt = 0;
int qh = 1;
while (qt < qh) {
int v = q[qt++];
for (int e = firstEdge[v]; e >= 0; e = nextEdge[e]) {
int u = edgeDst[e];
if (edgeFlow[e] < edgeCap[e] && level[u] < 0) {
level[u] = level[v] + 1;
q[qh++] = u;
}
}
}
return level[t] >= 0;
}

}

}

a.poorakhavan

Guest Blogger

categories & Tags

UNLEASH THE GIG ECONOMY. START A PROJECT OR TALK TO SALES
Close


text:
?

Free Version
Moderate

# Dueling Magnetic Fields, Direction

APPH12-VBWVK2

Two strong magnets are placed so that the magnetic field of magnet $A$ has a magnetic field of $0.80$ $T$ to the West at point $P$. Magnet $B$ has a magnetic field of $1.50$ $T$ to the North at the same point $P$.

What is the direction of net magnetic field at point $P$ due to both magnet $A$ and magnet $B$?

A

It cannot be represented by a single direction because magnetic fields do not interact.

B

32 degrees North of West.

C

28 degrees North of West.

D

62 degrees North of West.


text:
# Nusselt number

(Difference between revisions)
 Revision as of 08:37, 5 December 2005 (view source)← Older edit Revision as of 10:23, 5 December 2005 (view source)Jola (Talk | contribs) Newer edit → Line 3: Line 3: Nusselt number is dimensionless number which is useful for convective heat transfer calculations. heat transfer can be expressed in terms of Nusselt number. it can be defied as follows: Nusselt number is dimensionless number which is useful for convective heat transfer calculations. heat transfer can be expressed in terms of Nusselt number. it can be defied as follows: - '''[itex]Nu = h*l/k[/itex]''' + :[itex]Nu = h*l/k[/itex] Where Where

## Nusselt number

Nusselt number is dimensionless number which is useful for convective heat transfer calculations. heat transfer can be expressed in terms of Nusselt number. it can be defied as follows:

$Nu = h*l/k$

Where

$Nu = \mbox{Nusselt number}$
$h = \mbox{Convection heat-transfer coefficient (constant of proportionality ) (W/m2 K)}$
$l = \mbox{Characteristic dimension of th esolid object (m)}$
$k = \mbox{Thermal conductivity of the solid (W/m k)}$


text:
Browse Questions

# According to Einstein’s photoelectric equation, the plot of the KE of the emitted photoelectrons from a metal vs the frequency of the incident radiation gives a straight line whose slope

$\begin{array}{1 1} \text{(a) depends on the nature of the metal used} \\ \text{(b) depends on the intensity of the radiation} \\ \text{(c) both A & B} \\ \text{ (d) Depends on neither A nor B} \end{array}$

Ans : (d)
Einstein’s photoelectric equation is $KE_{max} = hc \: – \phi$……..(i)
The equation of the line is $y=mx+C$………..(ii)
comparing eqns. (i) & (ii),
$m = h, c = -\phi$
So, slope is a constant and independent of intensity of radiation
edited Mar 26, 2014